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Our recent quasi-two-dimensional thermodynamic description of thin-liquid films stabilized by col- 
loidal particles is generalized to describe nonuniform equilibrium states of films in external potentials 
and nonequilibrium transport processes produced in the film by gradients of thermodynamic forces. 
Using a Monte-Carlo simulation method, we have determined equilibrium equations of state for a 
film stabilized by a suspension of hard spheres. Employing a multipolar-expansion method combined 
with a flow- reflect ion technique, we have also evaluated the short-time film- viscosity coefficients and 
collective particle mobility. 



I. INTRODUCTION 

Owing to the presence of oscillatory structural forces, 
static and dynamic properties of thin liquid films stabi- 
lized by colloidal particles, micelles, or macromolecules 
are quite different than the properties of particle-free 
films P, [2) Q. A striking example of thin-film behav- 
ior that is caused by oscillatory structural forces is the 
stepwise-thinning (i.e., film stratification) phenomenon. 
Stepwise thinning — or thickening Q — occurs in liquid 
films stabilized by colloidal particles [1, d, B, Bl , poly- 
electrolytes [t| and surfactant solutions . Such 
films undergo a sequence of stepwise transitions between 
regions of different, but uniform, thickness. In some 
cases, quite complex multiphase structures are obtained 
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The stepwise-thinning in particle-stabilized thin li quid 
films has been studied experimentally d, 0, 0, [2 flUllrjl] 
and theoretically 0, fl HE E d, Hfl • Most of the anal- 
yses of the mechanism of stepwise thinning focused on 
the role of normal structural force produced by the sus- 
pended particles [!, d, [3 [H, [H, H3| • However, as we 
have pointed out in our recent paper [l7j . not only the 
normal but also lateral structural forces play an impor- 
tant role in the film-stratification phenomenon. 

To include the lateral structural forces in our descrip- 
tion of phase equilibria in particle-stabilized thin liquid 
films we have proposed a quasi-two-dimensional thermo- 
dynamic formalism (l7l |2l| . In this formalism the film 
area, film volume, and number of colloidal particles are 
independent parameters of state. The corresponding con- 
jugate intensive parameters are the film tension, normal 
pressure, and chemical potential of the colloidal particles 
in the film. 

The film tension involves the anisotropic part of the 
osmotic-pressure tensor of the suspension confined be- 
tween the film interfaces. Hence, the lateral component 
of the structural force produced by the suspended parti- 
cles contributes to mechanical equilibrium in the lateral 
direction (i.e., direction along the film). An analysis of 
the lateral structural forces is thus essential for describing 



film properties, in particular, for formulating the equilib- 
rium coexistence conditions between film phases of dif- 
ferent thickness [l7| . 

In the present paper the ideas introduced in [l7| are 
developed in two directions. First, we propose a theoreti- 
cal description of equilibrium states of films with varying 
thickness h and nonuniform particle number density per 
unit area n c when lateral external forces are present. Sec- 
ond, we present a systematic analysis of macroscopic film 
dynamics, based on effective 2D transport equations that 
involve viscous dissipation and diffusive particle flux. 

We also provide numerical results for equilibrium and 
nonequilibrium properties of a film stabilized by a sus- 
pension of hard spheres. The normal stress and film ten- 
sion are obtained from the contact value of the pair dis- 
tribution, which is calculated using a Monte-Carlo tech- 
nique. To determine the chemical potential of the par- 
ticles confined in the film we propose a new evaluation 
method, based on integration of the Gibbs-Duhcm re- 
lation that links variations of the particle chemical po- 
tential to the corresponding variations of the normal os- 
motic pressure and film tension. The transport coeffi- 
cients in our nonequilibrium theory are determined using 
a Stokesian-dynamics algorithm for a periodic system of 
spheres (23l . |24| . combined with a flow- reflect ion tech- 
nique. 

Our paper is organized as follows. A general frame- 
work for a description of nonuniform equilibrium states in 
particle-stabilized thin- liquid films is developed in Sec.lHl 
followed by an analysis of nonequilibrium transport pro- 
cesses in Sec. Mil Short-time transport coefficients are 
evaluated in Sec. IIV1 and concluding remarks are given 
in Sec. N\ 



II. EQUILIBRIUM DESCRIPTION 

A. Particle-stabilized films 

We consider a thin liquid film of thickness h and area 
A stabilized by a colloidal suspension of hard spheres 
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FIG. 1: Geometry of a particle-stabilized thin liquid film. 



of diameter d ~ h. The film interfaces are surfactant- 
free and have a constant interfacial tension a. The film 
is surrounded by an inert gas of a constant pressure 
p^ . The viscosity and mass density of the gas are much 
smaller than the corresponding parameters for the fluid 
in the film; therefore, the gas interacts with the film only 
through the static pressure p^ . The temperature T in 
the system is uniform. 

We assume that the film is approximately planar. The 
film thickness and the number of particles per unit area 
are either constant or slowly varying in the lateral di- 
rections x and y on the length scale L 3> h. Under 
these conditions the transverse suspension structure re- 
laxes on the time scale t± that is much shorter than the 
time scale £|| for evolution of the long-wavelength lateral 
modes. After the time tj_, the suspension in the film is 
thus in a local-equilibrium state, characterized by den- 
sities per unit area of conserved quantities such as the 
volume or number of colloidal particles. 

Three-dimensional macroscopic quantities that appear 
in our analysis (e.g., stress and pressure tensors, suspen- 
sion velocity, and particle flux) are the ensemble averages 
of the corresponding microscopic quantities. Assuming 
that the film properties are slowly varying in the lateral 
directions, the local area averages and ensemble averages 
are equivalent. In the effective quasi-two-dimensional de- 
scription of the film, the averages across the film of 3D 
quantities are evaluated, in addition to the ensemble av- 
erages. 

In all our numerical examples the suspending fluid is 
treated as an incompressible continuous medium. The 
fluid affects the thermodynamic state of the system only 
through the isotropic contribution to the pressure tensor 
and through the constant- volume constraint. The par- 
ticles in the film interact via the hard-sphere potential. 
We note, however, that our thermodynamic considera- 
tions apply more generally, e.g., to suspensions of charged 
particles or small particles with finite-size hydration lay- 
ers. In our hydrodynamical calculations the creeping-flow 
conditions are assumed, and many-body hydrodynamic 
interactions of the particles in the film are fully taken 
into account. 



B. Normal pressure and film tension 

Unlike the corresponding quantity in bulk suspensions, 
the equilibrium pressure tenor p in particle-stabilized 
films is anisotropic. This anisotropy stems from different 
particle ordering in the transverse and lateral directions. 
For films with parallel (or nearly parallel) planar inter- 
faces, the pressure tensor has only the normal and lateral 
components, 

P = P±e z e z +p||l s , (1) 

owing to the axial symmetry of the problem. Here l s = 
&x&x + B y e y is the lateral unit tensor (where and ej, 
are unit vectors in the lateral directions x and y) , and e 2 
is the unit vector in the transverse direction z (normal 
to the film). 

Due to continuity of the momentum flux, the normal 
pressure p± is independent of the transverse coordinate 
z. However, the lateral component pii of the pressure 
tensor depends on the transverse position in the film. In 
our paper the average value across the film of a quantity 
b (such as the lateral stress pn) will be denoted by the 
overbar, i.e., 

b = h~ 1 / dzb(z). (2) 
Jo 

Normal and lateral stresses acting on a section of the film 
are schematically illustrated in Fig. [TJ 

The infinitesimal work associated with the change of 
the film thickness h and area A is 

$W = -Ap± Ah + (2cr - hp\\ ) dA. (3) 

This expression includes the pressure and surface-tension 
contributions. 

The quantity /|| = (2cr — hp\\ ) in Eq. ([3]) represents the 
lateral force (per unit length of a film section) resulting 
from the interfacial-tension and pressure contributions. 
This force, however, cannot be interpreted as film ten- 
sion, and a correct relation for the film tension is given in 
Eq. (J5J. The term (2cr - hp\\)dA in Eq. © describes the 
work associated with a change of the film area at a con- 
stant thickness h. The description in terms of indepen- 
dent variables h and A would give meaningless results for 
a film composed of an incompressible suspension, because 
any area change in such a system must be accompanied 
by the corresponding thickness change. Moreover, since 
we have pii = p± for a film composed of a simple isotropic 
fluid, the force fu depends on the external pressure p^ 
(which is in equilibrium with the normal pressure p± in 
the film). The film tension, however, should be pressure 
independent. 

A correct film-tension representation of the equilibrium 
states of the film is obtained by rewriting the pressure 
tensor (JTJ in the form 

p =p±\ + (p\\ -p±)\ s , (4) 



3 



where I = \ s +e z e z is the 3D unit tensor. The first term on 
the right-hand side of the above relation is an isotropic 
pressure component that should not contribute to the 
film tension (as discussed above). The second term, inte- 
grated across the gap, corresponds to the excess lateral 
force per unit length of film section. Only this excess 
force contributes to the film tension 7. Combining the 
excess lateral pressure with the surface-tension contribu- 
tion we obtain the relation for the film tension, 

7 = h{px -P\\) + 2er (5) 

and the corresponding expression for the work, 

aw = ~ P± dV + jdA, (6) 

where V = Ah is the film volume. The term 7 a A in 
the above equation describes the work associated with 
the change of the film area at a constant film volume, in 
contrast to the corresponding term in Eq. (|3"|). 

In mechanical equilibrium, the normal pressure in the 
film equals the constant pressure of the surrounding gas, 

P±=P (C) - (7) 

The excess lateral stresses have to balance as well, which 
implies that the film tension is position independent 

7 = const, (8) 

provided that there are no external forces acting on the 
film. 

In the above discussion we assumed that there is no 
direct interaction between the film interfaces. However, if 
there is such an interaction force (e.g., the van der Waals 
attraction or electrostatic repulsion), the normal pressure 
p± in Eqs. I©-® should be replaced with p' ± = p± + f n , 
where f n is the normal force per unit area, and f n < 
corresponds to attraction. 



C. Quasi- two-dimensional thermodynamic 
description 

Owing to the time scale separation t± <C tu , the relax- 
ation of the suspension structure in the transverse direc- 
tion is much faster than the relaxation of long- wavelength 
lateral modes. In a quasi-two-dimensional thermody- 
namic description [l7| all details of the transverse struc- 
ture of the film are thus averaged out. The film thickness 
h = V/A, however, is retained as a thermodynamic vari- 
able, because it can be controlled by varying the external 
pressure. 

The fundamental thermodynamic relation for the free 
energy F in the film-tension representation can be ob- 
tained from the expression for work ([6]) by using the 
standard entropy-maximization arguments [251 ]. Treat- 
ing the suspension in the film as a two-component fluid 
(with the colloidal particles regarded as macromolecules) 



we get the fundamental relation in the free-energy repre- 
sentation 

dF= -S'dT-p_ L d^ + 7dA + ^ c diV c + ^ f d7V f , (9) 

where S is the entropy, and N c and Nf denote, re- 
spectively, the number of colloidal particles and solvent 
molecules in the system, and /i c and fj,{ are the corre- 
sponding chemical potentials. According to the above re- 
lation and the free-energy-minimum principle, a particle- 
stabilized film in thermodynamic equilibrium must sat- 
isfy the thermal equilibrium condition T = const, the 
mechanical equilibrium conditions ([7]) and l[5]). and the 
chemical equilibrium condition 

/i c = const, /if = const. (10a, b) 

A direct consequence of the fundamental relation © 
is the Gibbs phase rule for coexisting regions of different 
thickness h in stratified films in thermodynamic equilib- 
rium. Assuming that the suspending fluid behaves as 
incompressible continuum medium, its degrees of free- 
dom can be neglected (as discussed in Sec. Ill El) . Thus, 
the state of each phase is described by three densities 
per unit area: the excess entropy density, the film thick- 
ness h = V/A and the number density of the colloidal 
particles n c — N c /A. There are also four equilibrium 
conditions corresponding to the intensive parameters T, 
p±, 7, and p c . It follows that a two-phase system has 
/ = 2 degrees of freedom, and the maximum number of 
coexisting phases is k = 4. 

We note that in experiments with vertical films, up 
to seven coexisting phases of different thickness were ob- 
served 0. In such systems, additional thermodynamic 
degrees of freedom are provided by the external grav- 
ity potential, and it is possible that not all equilibrium 
conditions are satisfied (due to dynamical constraints re- 
sulting in slow relaxation of certain intensive-parameter 
differences [l7|). 



D. Nonuniform systems 

1. External forces 

The physical meaning of film tension ([5]) can be further 
elucidated by considering nonuniform equilibrium states 
of a film in external lateral force fields. In Sec. IIIII an 
analysis of such nonuniform film states will also help to 
determine the linear constitutive relations governing film 
dynamics. 

We consider potential forces 

F f = -V„* f> F c = -V s * c , (lla,b) 

acting, respectively, on the fluid molecules and colloidal 
particles. It is assumed that the potentials 'ff and \& c 
depend only on the lateral position p — xe x + ye y . 
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The force fields (fTTjl are normalized per single particle. 
Thus, the total force density per unit area of the film is 

f = n f F f + n c F c , (12) 

where Uf and n c are the local values of the area number 
densities fif — Nf/A and n c — N c /A of the fluid molecules 
and colloidal particles at the lateral position p. Similarly, 
the total force density per unit volume is given by the 
relation 

f = ?i f Ff + n c F c , (13) 

where rif and n c are the local number densities per unit 
volume of the fluid molecules and colloidal particles at a 
(three-dimensional) position r. 

2. Lateral stress balance 

The balance between the film tension and external 
force (fT2]) can be obtained from the full 3D stress-balance 
equation 

V-p = f (14) 

(which is valid in the whole space). In the region inside 
the film, the pressure tensor p is given by Eq. ([4|) ; outside 
the film it is equal to the external pressure tensor, 

p=p (c) l. (15) 

In Eq. IHD , the interfacial contribution is neglected un- 
der the assumptiona that there is no spacial variation 
of the interfacial tension, a = const, and that the film 
thickness varies on the length scale L ^> h (and thus the 
film curvature is negligible). 

Combining relations (|4]) and (fl~5|) for the pressure ten- 
sor inside and outside the film, and using the normal- 
stress boundary condition ((TJ, we find that the divergence 
of the isotropic pressure component vanishes in the whole 
space. The integral of the remaining part with respect 
to the variable z yields the lateral stress balance (cf. Ap- 
pendix [K§ 

V s -h(p ]l -px)\ s = f, (16) 

where f is the area force density (|12j) . By using definition 
([5]) of film tension, and applying the assumptions L 3> h 
and a = const to include the interfacial-tcnsion term, the 
lateral stress-balance equation (fT6|) can be expressed as 
the hydrostatic condition for the film-tension gradient 

V s7 = -f. (17) 

Equation (|17[) requires that, in equilibrium, the area 
force density be a gradient of a 2D potential, 

f = -V B f (18) 



Thus, in the presence of external force field (fl"2f . the 
equilibrium condition (jHJ) is replaced by the hydrostatic 
condition for the film tension, 

7 — -ip = const. (19) 

A more general discussion of the lateral stress balance 
in the film, applicable to equilibrium and nonequilibrium 
states, is given in Appendix [Al 

3. Chemical potentials 

External forces (fTTj) modify not only the lateral stress 
balance, but also the equilibrium conditions (TIT)]) for the 
chemical potentials of the suspending fluid and colloidal 
component. As in the standard thermodynamics of bulk 
systems, the modified chemical equilibrium conditions 
are 

/if + \&f = const, fi c + ty c = const. (20a, b) 

We note that the above expressions are consistent with 
the lateral stress balance (fT7| . owing to the Gibbs- 
Duhem relation 

- SdT + Vd Pl _ - Ad-f-N c dfi c - N f dfM = (21) 

(which is obtained from @ by the usual Legendre trans- 
formation of all extensive variables). At constant T and 
px, the Gibbs-Duhem relation yields 

V s 7 = -n c V s /i c - UfV s /if. (22) 

By inserting the equilibrium conditions (|20p into (|22|) and 
using the definition (| 1 2[) of the area force density f , the 
lateral stress-balance equation (fT7|) is recovered. 

E. Incompressible suspension 

In our further considerations we focus on films com- 
posed of an incompressible suspension. It is assumed 
that both the solvent and the suspended particles are in- 
compressible, and changes of the suspension volume due 
to interaction of particle hydration layers are negligible. 

The number of solvent molecules in such a film is set by 
the film volume V and the colloidal-particle volume frac- 
tion 4> c — ^nd 3 fi c , where n c = N c /V is the particle num- 
ber density averaged over film thickness. The chemical 
potential fi{ of the solvent is not an independent function 
of state, either. Under isothermal conditions, assumed 
herein, /if is fully determined by the Gibbs-Duhem rela- 
tion for the suspending fluid 

n[ 0) d/if = dp t (23) 

(where n[ ^ is the number density per unit volume for 
the fluid molecules in pure solvent). The thermodynamic 
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degrees of freedom associated with the solvent variables 
/Xf and Nf can thus be eliminated from the fundamental 
relation 0. 

The thermodynamic description of an incompressible 
system can be rephrased in terms of appropriate excess 
quantities describing the colloidal contributions to ther- 
modynamic functions of state. In this section we define 
excess quantities and formulate the corresponding equi- 
librium conditions. 

The pressure tensor in the film, 

P=Pfl + Pc (24) 

is decomposed into the isotropic fluid pressure pi and the 
osmotic pressure tensor 

p c = p c ±e z e z +_p c ||l s (25) 

that results from the presence of the colloidal particles. 
The isotropic fluid pressure pt does not contribute to film 
tension, according to the definition ([5]). It is, therefore, 
convenient to rewrite the film tension 

7 = 7c + 2a (26) 

as a sum of the interfacial tension term and the colloidal 
contribution 

7c = h(p c ± -p c \\). (27) 

Since we assume no variation of the interfacial tension 
along the film, the colloidal contribution 7 C differs from 
the full film tension 7 only by a constant. Thus, the 
hydrostatic-equilibrium condition for film tension (|19[) re- 
mains, essentially, unchanged, 

7 C — 1/; = const. (28) 

In a system with no external potentials, the fluid pres- 
sure satisfies the equilibrium condition pf — const. In 
the presence of an external potential (Tib), the chemical- 
equilibrium condition (|20h ) and the Gibbs-Duhcm rela- 
tion ([23D yield 

p { + fa = const, (29) 

where 

^ { = nf ) ^ l (30) 

is the potential of the volume force acting on the fluid 
phase. Equation (j2"9"| represents the usual hydrostatic- 
equilibrium condition for the suspending fluid. 

With the help of the hydrostatic relation (f29|) . the 
constant-normal-pressure condition ([7]) can be trans- 
formed into the corresponding equilibrium condition for 
the osmotic normal pressure p c ±, 

p c ± — ipf = const. (31) 

The equilibrium stress-balance equation (|3"Tj) relates the 
normal osmotic pressure in the film to the potential ipt of 



the lateral force acting on the fluid. Thus, in equilibrium 
states with a nonzero lateral force, the normal osmotic 
pressure varies along the film. 

The remaining equilibrium condition that has to be 
considered, is relation l|20b) for the chemical potential 
of the colloidal particles /i c . In the fundamental relation 
© , the term /z c dJV c corresponds to the change of the 
free energy F when a colloidal particle is added to the 
system at a fixed volume and number of fluid particles. 
However, in an incompressible suspension, a particle that 
is moved from one subsystem of the film to another is al- 
ways replaced by the equivalent amount of the suspend- 
ing fluid. Accordingly, the quantity /x c involves not only 
the configurational contribution resulting from the dis- 
tribution of the colloidal particles, but also an additional 
part corresponding to the amount of the fluid displaced 
by a particle. 

To define the excess chemical potential p,' c that is free of 
the fluid contribution, we begin with the Gibbs-Duhcm 
relation (|2"Tj) . By combining (|2"Tj) with the Gibbs-Duhcm 
relation for the suspending fluid (|2"3")l we obtain the iden- 
tity 

Vd Pl _ - Ad 1 - N c d(i c - (1 - <t> c )Vd Pi = 0, (32) 

where the constant-temperature conditions are assumed. 
Defining the excess chemical potential of the colloidal 
particles as 

/4 = Mc - VpPf (33) 

(where v p is the volume of a colloidal particle) , and using 
relation pj_ = Pc ±+Pf, we find the Gibbs-Duhem relation 
for /4, 

dfj,' c = h^ 1 dp cl ^ - h^ 1 d7 c . (34) 

The above equation indicates that the quantity /j,' c is de- 
termined (up to a temperature-dependent additive con- 
stant) by the osmotic normal pressure and film tension, 
which justifies our interpretation of /i c as the excess 
chemical potential. 

In order to rephrase the chemical equilibrium condition 
for the colloidal particles in terms of the excess chemical 
potential p,' c , relations (|2^)) and (|20b ) are inserted into 
(|3"3"1) . As the result we get the chemical equilibrium con- 
dition for the excess quantities 

n' c + * c = const, (35) 

where the excess potential 

% = * c - Vpt/H (36) 

of the external force acting on the colloidal particles in- 
cludes the buoyancy contribution resulting from the pres- 
sure gradient in the suspending fluid. 
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F. Numerical results 

1. Evaluation technique 

In this section we present some numerical results for 
the osmotic-pressure tensor p c , film tension j c , and parti- 
cle chemical potential /i' c for a film stabilized by a suspen- 
sion of hard spheres of diameter d. A standard Metropo- 
lis Monte-Carlo algorithm was used to generate equilib- 
rium ensemble for a system of a hundred hard spheres in 
a film with periodic boundary conditions in the lateral 
directions. The results were obtained for several volume 
fractions of colloidal particles <f> c , as a function of film 
thickness h. 

a. Osmotic-pressure tensor The ensemble-averaged 
osmotic pressure tensor has been evaluated using a hard- 
sphere generalization of the Kirkwood-Buff expression 
for the pressure tensor in an inhomogeneous fluid [2t| . 
Accordingly, the colloidal contribution to the pressure 
tensor p c , averaged over the space V = A(h— d) accessible 
to the particle centers, is obtained from the relation 




Here n^ri,^) is the two-particle equilibrium reduced 
distribution, f = r/r is the unit vector along the line 
passing through the particle centers, the integration is 
over the contact configurations, (■■■)■& denotes the vol- 
ume average over the region V, and n c — N c /V is the 
corresponding particle number density. 

Equation (|37[) generalizes a well-known expression for 
the pressure in a bulk hard-sphere system in terms of 
the contact value of the radial distribution. The contact- 
value relation (|37|) can be derived either from the col- 
lisional contribution to the momentum flux in a hard- 
sphere system, or by passing to the hard-sphere limit in 
the Kirkwood-Buff expression [26[ for the pressure tensor 
in an inhomogeneous fluid. 

In our numerical calculations, the average (|3T[) is eval- 
uated from a sum of contributions originating from pairs 
of particles separated by a distance < d(l + e), where 
e <C 1. The results for finite values of e are extrapolated 
to obtain the limiting result for e — > 0. 

The contact-value relation (]3"T]) provides a convenient 
basis for numerical evaluation of the pressure tensor in 
the film. However, p c involves only averaging over the 
accessible volume V , and does not incorporate the mo- 
mentum transfer from the particle centers to the wall. 
To obtain the pressure tensor p c averaged over the whole 
volume of the film V, we need to consider this momentum 
transfer. 

The normal momentum transfer occurs when a par- 
ticle collides with the interface [271 ] . Thus, there is a 
constant normal momentum flux in the geometrically ex- 
cluded regions of thickness a = ^d next to each interface. 
It follows that the normal component of the pressure is 
not affected by the choice of the region (V or V) over 



which the average is performed. In contrast there is no 
lateral momentum transport in the excluded layers next 
to each wall, and the averages over the volumes V and V 
differ by the factor V/V = (h- d)/h. 

The above arguments imply that the normal and lat- 
eral components of the tensor p c averaged over the vol- 
ume V are 

Pc±=Pc±, p c \\ = /i _1 - d)p c \\, (38) 

where p c = p c ±e z e z + p c u\ s is given by the contact-value 
relation ({57]) . The components (|38[) are used to determine 
the colloidal contribution to film tension (]27[) . 

b. Chemical potential Our quasi-two-dimensional 
thermodynamic formalism can be conveniently used for 
evaluation of the chemical potential of the colloidal parti- 
cles in the film [i' c from the results for the normal pressure 
and film tension. Our method relies on numerical inte- 
gration of the Gibbs-Duhem relation (]34]) , using results 
of Monte-Carlo simulations at different values of particle 
volume fraction <j) c and film thickness h. 

For convenience, in our calculations we use a modi- 
fied Gibbs-Duhem relation based directly on the pressure 
tensor ([57]) . without introducing the averages (|3"5)) over 
the full film volume V. In this modified formulation, the 
Gibbs-Duhem relation (]34]) is replaced with 

d/4 = n^ 1 dp cl _ - h^ 1 d%, (39) 

where the modified film tension is given by the formula 

7c = (ft-d)(Pe±-Pc||) > (40) 

instead of (]27j) . 

We note that equations ([34]) and |39]) are equivalent, 
which follows from the identities 

7c = 7c + dp c± , (41) 

n-^n-i + dh- 1 , (42) 

where, as previously, the overbar and the hat indicate 
the average particle density in the volumes V and V, 
respectively, and the tilde indicates the particle density 
per unit area. 

We have verified that numerical integration of Gibbs- 
Duhem relation (|39p over different paths in the thermo- 
dynamic space of states yields equivalent results. To ob- 
tain proper normalization of the chemical potential, rela- 
tion (|39p is integrated from a low-density state, where we 
can use the perfect-gas result for the chemical potential. 

Since equilibrium properties of a confined suspension 
do not depend on the details of the dynamic boundary 
conditions at the interfaces of the film (e.g., free interface 
or rigid wall), our method for evaluating the chemical po- 
tential can be applied not only to thin liquid films, but 
also to other systems with a slit-pore geometry. In par- 
ticular, our method can be used to determine the chem- 
ical potential of spherical macromolecules in the gap be- 
tween two much larger particles in a system where de- 
pletion interactions are important. Our method can also 
be generalized to systems interacting via continuous or 
anisotropic potentials. 
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FIG. 2: (a) Normal and (b) lateral osmotic pressure in the 
film, versus film thickness h, for particle volume fractions 4> c = 
0.25 (dotted line), 0.3 (dashed), 0.35 (dot-dashed), and 0.4 
(solid). 



2. Results 
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FIG. 3: Particle contribution to film tension, versus film 
thickness h, for the same values of volume fraction <j> c as in 
Fig.H 



leading to the pressure oscillations and of the source of 
the shift between the maxima of p c ± and pd| can be 
gained from an analysis of the suspension structure for 
different values of the film thickness h. In Fig. [5] we show 
the local particle density n c in the film for two values 
of interface separation, corresponding to the maximum 
(h/d = 1.6) and minimum (h/d = 2.25) of the normal 
pressure p c ± [cf., Fig. H£a)]. 

The results shown in Fig. [5] indicate that for both val- 
ues of h/d, there are two particle layers in the film. For 
h/d = 1.6, the two layers barely fit into the available 
space V; particle concentration is thus strongly focused 
in thin regions adjacent to each interface, and the sus- 
pension exerts high normal pressure p c ±_ on the interfaces 
of the film. For h/d — 2.25, the particle layers are more 



The dependence of the normal and lateral pressure 
components (|38|) on the film thickness is illustrated in 
Fig. [21 which shows p c ± and p c \\ versus h for several val- 
ues of the particle volume fraction cf> c . The corresponding 
behavior of the film tension j c and excess chemical po- 
tential fi' c is depicted in Figs. l3l and l4l [28j]. 

The results plotted in Fig.[2]indicate that bothp c ± and 
p c M are oscillatory functions of h. However, the phase of 
the oscillations of the normal and lateral pressure com- 
ponents significantly differ, with the maxima of p c u oc- 
curring at smaller values of h than the maxima of p c ±. 
Film tension (|27p and chemical potential fi' c [which is ob- 
tained from the pressure tensor by integrating the Gibbs- 
Duhem relation (f39j) ] are also oscillatory functions of h. 
(cf., Figs. [3] and HJ). 

The oscillatory behavior of the pressure tensor is typ- 
ical of effective structural interactions in confined do- 
mains. The nonmonotonic dependence on h stems from 
the short-range layered ordering of the suspension in the 
film [29j |. A qualitative understanding of the mechanism 







/ /N \ 
















• , 







h/d 

FIG. 4: Excess chemical potential of the particles in the 
film, versus film thickness h, for the same values of volume 
fraction <f> c as in Fig. [2] The normalization corresponds to the 
asymptotic behavior fi' c = kTln(d 3 h c ) for cf> c <JC 1. 
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FIG. 5: Local particle number density n c (per unit volume), 
versus position z across the film, for volume fraction cf> c = 0.4 
and two values of film thickness h/d = 1.6 (solid line) and 
h/d = 2.25 (dashed line). Relative to the film thickness, the 
particle distribution is much narrower for h/d = 1.6 (maxi- 
mum of p c _L in Fig. [2} than the distribution for h/d — 2.25 
(minimum of p c ±)- 



diffuse, and the normal pressure is much lower. 

In contrast to the behavior of p c ±, the lateral pressure 
component p c || remains moderate both for h/d = 1.6 and 
h/d = 2.25 [cf. the results shown in Fig. [2{&)]. This is 
because for h/d = 1.6 there are very few particles in the 
central portion of the film, which implies that there is 
no lateral momentum transport in this region. The av- 
erage value thus remains relatively low, in spite of high 
momentum flux in the layers of high particle concentra- 
tion. As a result of this mechanism, the maxima of the 
lateral pressure p c \\ (averaged across the film) occur for 
interface separations where both, the local lateral mo- 
mentum transport is relatively large, and it takes place 
over sufficiently broad regions. 

The oscillatory form of the intensive thermodynamic 
parameters p c ±, 7 C , and p' c strongly affects properties 
of particle-stabilized films. First, owing to these oscil- 
lations, the coexistence between uniform film phases of 
different thickness is possible, as observed experimentally 
d, H, 0, HI and discussed in our recent paper [l^] ■ 

Next, a small perturbation 6h of the film thickness 
from a uniform stable equilibrium state produces restor- 
ing force proportional to the lateral gradient of Sh. In 
contrast, in particle-free films the leading-order restor- 
ing force results from the capillary pressure, and it is 
thus proportional to the gradient of the film curvature 
k ~ Vg/i. For long-wavelength perturbations, the restor- 
ing force acting in particle-stabilized films is thus sig- 
nificantly stronger than the restoring force in particle- 
free films. Therefore, particle-stabilized films often re- 
tain nearly uniform thickness during the drainage pro- 
cess (except for stepwise changes associated with phase 
transitions), whereas particle- free films usually assume a 



dimpled shape with nearly constant curvature. 

III. MACROSCOPIC TRANSPORT PROCESSES 

After discussing equilibrium thermodynamics of 
particle-stabilized thin liquid films, we now focus on lin- 
ear transport processes that occur when such films are 
perturbed from equilibrium. In Sec. IIII Al we present 
quasi-two-dimensional continuity equations for conserved 
parameters characterizing local-equilibrium states of the 
film, and in Sees. IIII Bl and IIII CI we formulate constitu- 
tive relations between thermodynamic forces and macro- 
scopic fluxes in the conservation equations. 

A. Conservation equations 

The thermodynamic analysis presented in Sec. |TT] indi- 
cates that a particle-stabilized thin liquid film behaves as 
a two-component 2D compressible fluid, even though the 
suspension in the film is incompressible. The film thick- 
ness plays the role of the mass density, per unit area, 
and the number of colloidal particles per unit area corre- 
sponds to the concentration of the solute. 

In the quasi-two-dimensional description, the macro- 
scopic motion of the suspension in the film is described 
by the lateral velocity 

v s = l s -v, (43) 

where v denotes the ensemble-averaged velocity field, and 
the overbar denotes the average across the film ©. The 
macroscopic motion of the particles in the film is charac- 
terized by the 2D lateral particle flux 

j c = l 8 -/ijc, (44) 

where j c is the 3D particle flux j c averaged across the 
film. We use the normalization where the particle flux j c 
is measured with respect to the local suspension velocity 
v s . 

Definition (|4"3"1) and the incompressibility condition 
V • v = for the 3D flow field yield the 2D continuity 
equation for the film thickness h, 

dh 

— = -V s .(hv s ), (45) 

where t denotes time. Similarly, the evolution of the 
particle density per unit area is governed by the 2D con- 
tinuity equation 

^ = -V s -(n c v s +j c ). (46) 

To obtain a full description of the film dynamics, the 
continuity relations (|45[) and (|46|) have to be supple- 
mented with the stress-balance equations and the con- 
stitutive equations for the stress tensor and the particle 
flux in the film. 
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B. Stress-balance equations 

1. Lateral stress balance 

According to the results given in Appendix [A] the 
stress-balance equation in the film-tension representation 
involves the effective 2D surface stress f and area-force 
density (O, 

V s -f = -f. (47) 

The effective surface stress, defined by relations (|A7[) 
and (|A13[) . equals the excess stress beyond the constant 
isotropic background pressure p^\, integrated across the 
film. As in all our dynamic considerations, in the stress- 
balance equation (|47[) . creeping- flow conditions are as- 
sumed. 

In equilibrium, the 2D surface-stress tensor reduces to 
the isotropic film tension, r = j\ s . In a nonequilibrium 
state, the surface stress 

t = 7 I. + f ' (48) 

is a sum of the local-equilibrium contribution 7l s and the 
viscous part f . In the above equation, the local value 
of the film tension 7 is given by the equilibrium equation 
of state, consistent with the standard local-equilibrium 
assumption [30| . For a particle-stabilized film, the film 
tension 7 = j(h, n c ) is fully determined by the local film 
thickness h and local number of colloidal particles per 
unit area of the film n c . 

In our further discussion, we assume the linear- 
response regime and long- wavelength limit. Under these 
conditions, the viscous stress f' is related to the lateral 
velocity (|4"3")l via two film-viscosity coefficients: the shear 
and the expansion viscosities fj and it. Accordingly, we 
have 

t' = 2r)[V s v s ] d + K(V s -v s )l s , (49) 

where 

[V s v s ] d = i[V s v s + (V s v 8 )t - V s • v s |J (50) 

is the deviatoric part of the lateral strain-rate tensor 
V s v s , with the dagger denoting the transpose. 

2. Normal stress balance 

Owing to the suspension incompressibility, the normal 
balance of the total stress in the film does not have an in- 
dependent dynamical meaning, unlike the lateral balance 
(|47|). However, the normal osmotic pressure component 
is an important dynamic quantity, because it drives a 
lateral particle flux. These results are derived below. 

In a nonequilibrium state of the film, the equilibrium 
condition ([7|) is replaced by the normal-stress balance 

Ti = - P (G) , (51) 



where t± is the normal component of the stress tensor r. 
Unlike the effective surface stress (|48|) - the normal stress 
t± cannot be uniquely decomposed into the equilibrium 
and nonequilibrium components. This is because in an 
incompressible suspension the isotropic component of the 
pressure tensor (and thus also the normal pressure p±) 
is not a well-defined function of the densities of the con- 
served quantities h = V/A and n c = N c /A. There is thus 
no independent constitutive relation for the nonequilib- 
rium part of tj_, and any unbalanced normal stress is 
immediately compensated for by a proper redistribution 
of the isotropic component of the 3D pressure tensor. We 
note that such redistribution does not modify the lateral 
stress balance in the film, because the surface stress f 
involves only differences between the normal and lateral 
components of the 3D stress, according to Eq. (|A13|) . 

While the normal stress cannot be uniquely decom- 
posed into the equilibrium and nonequilibrium parts, 
there exists another useful decomposition. Namely, the 
stress t j_ can be decomposed into the equilibrium osmotic 
normal pressure p c ± and the remaining nonequilibrium 
part t", 

T_L = Pc± + t". (52) 

The nonequilibrium stress component r" involves the 
3D isotropic fluid-pressure, and anisotropic contributions 
that result from (a) the viscous fluid flow and (6) devia- 
tion of the particle distribution from equilibrium. Owing 
to the equilibrium condition (|5ip , for given p^ the stress 
t" is entirely determined by the osmotic normal pressure 
PcA_, which is a well-defined function of h and n c . 

C. Particle flux 

According to the principles of the thermodynamics of 
irreversible processes [30| . in the linear-response regime 
the diffusive flux of colloidal particles is driven by gradi- 
ents of scalar intensive thermodynamic parameters that 
characterize the thermodynamic equilibrium conditions. 
In a particle-stabilized thin liquid film, the relevant in- 
tensive parameters are the normal osmotic pressure p c ± 
(as discussed above) and the chemical potential of the 
suspended particles fi' c . If a particle-stabilized film is 
subject to external lateral forces, we also have thermo- 
dynamic forces associated with the external potentials ipf 
and ^' c , according to the equilibrium conditions (|3ip and 
([55)) . However, the spatial variations of the film tension 
7 (playing a role analogous to the isotropic equilibrium 
pressure in simple fluid) produce only the macroscopic 
motion of the film. Hence, in the linear-response regime, 
we assume the following constitutive relation for the 2D 
diffusive particle flux 

Ic = -A„V s (/4 + K) - A P V s (p c± - Vf). (53) 

Using the Gibbs-Duhem relation (|34|) it can be shown 
that the kinetic coefficients A M and A p cannot be deter- 
mined separately. Therefore, we combine them into a 
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single mobility coefficient v' . We first consider a quies- 
cent film with no viscous stress (|4"9")> . In this case the 
variations of the osmotic contribution to film tension -f c 
are related to the external potentials #' c and i/'f via the 
hydrostatic condition 



V s7c 



which follows from Eqs. {28]), d30j) , and (36]). 
Gibbs-Duhem relation (IMl) we also have 



(54) 



From the 



V S 7 C = n c V s /4 - hV s p c ±. 



(55) 



Combining Eqs. ([53"]) and (|54"]) yields one-to-one corre- 
spondence between the thermodynamic forces on the 
right-hand side of the constitutive relation (fS"3|) . 



n c V s (/4 + ^) = V s (p c 



(56) 



Using (|56p to eliminate the osmotic normal pressure p c ± 
and potential ipt from Eq. (|53|) . we find 



where 



X^ + h 1 X V 



(57) 



(58) 



is the collective mobility of the colloidal particles in the 
film. We note that Eq. ([5T|) relates the particle flux j c 
to the excess quantities ((33)) and (|36|) . indicated by the 
prime in our notation for the mobility coefficient (|58|) (cf. 
Sec. HVA3l) . 

In noncquilibrium states with nonvanishing gradients 
of the film velocity v s , the balance equation for the lateral 
stress involves the viscous contribution V s • r , accord- 



ing to equations (|4"T]) and (|4"8|) . As a result, a modified 
relation (|56p between the thermodynamic forces in Eq. 
(|53[) includes a V s • t correction. This correction, how- 
ever, is small in the long-wavelength regime, because it 
depends on second derivatives of the film velocity v s , ac- 
cording to the constitutive relation for the stress (|49|) . It 
follows that in the linear-response and long-wavelength 
limit, particle transport in the film can be described us- 
ing the constitutive equation (|57[) with a single mobility 
coefficient v' . 



short-time transport coefficients also describe the linear 
response to high-frequency oscillatory forcing. 

In the short-time (or high-frequency) regime the char- 
acteristic timescale of the forcing is much shorter than 
the structural relaxation time of the suspension in the 
film. Thus the deviation of the particle distribution from 
equilibrium is insignificant, and the Brownian contribu- 
tion to the transport coefficients can be neglected. Ac- 
cordingly, evaluation of the transport coefficients requires 
solving appropriate many-particle hydrodynamic prob- 
lems for the equilibrium particle distribution of spheres 
in the film under creeping- flow conditions. 



1. Shear viscosity of the film 

To evaluate the effective shear viscosity coefficient fj, 
we consider a film undergoing the planar straining flow 
with the strain-rate tensor 



E= i[Vv+(Vv) t ] 



given by 



10 
0-1 




(59) 



(60) 



where the coordinates x and y are in the directions along 
the film. By the axial symmetry and translational invari- 
ance with respect to lateral directions, the corresponding 
ensemble-averaged 3D viscous-stress tensor has the form 



2t?ii(z)E, 



(61) 



where the mesoscopic viscosity coefficient rj\\ (z) depends 
on the transverse position in the film z. Since equation 
(f6"T) involves only the lateral components of the stress, 
the effective shear viscosity fj is obtained by averaging 
Eq. ([6"T|) across the film. Using notation ([2|) we thus get 



fj = hr)\\, 

which follows from expressions 
(TAT31) . 



35b, EI 



(62) 
7l), and 



2. Expansion viscosity of the film 



IV. EVALUATION OF THE SHORT-TIME 
TRANSPORT COEFFICIENTS 

A. Elementary hydrodynamic processes 

In this section we evaluate the effective film-viscosity 
coefficients fj and k and the collective particle mobility. 
The film interfaces are assumed to be surfactant free, so 
free-interface boundary conditions are assumed. 

We focus on the short-time transport coefficients cor- 
responding to the initial system response to a startup of 
a flow or a sudden application of an external force. The 



The effective expansion viscosity re is obtained by eval- 
uating the stress response of the film to the axisymmetric 
straining flow with the strain-rate tensor given by 



1 
1 




(63) 



which corresponds to a uniform expansion of the film. By 
axial symmetry, the 3D viscous-stress tensor can be ex- 
pressed as a combination of the normal and lateral terms, 



Kj_e z e 2 



c l|(^ 



(64) 
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and E is the strain-rate tensor (|63[) . Comparing relations 
([63)) and (|6T|) , we find that the expansion film viscosity 
and the amplitude (3 of the deviatoric part of the average 
stress (frjfjj) are related by equation 



R = 3hf3. 



3. Particle mobility 



(68) 
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FIG. 6: Reflection method: (a) single free interface; (6) rep- 
resentation of a particle-stabilized film in terms of an un- 
bounded periodic system. 



The collective mobility coefficient ([58)) was determined 
by evaluating the particle flux j c produced by a constant 
lateral external force (jllb ) in a uniform film with zero net 
force ([T2"[) and zero net volume flow. The zero-net-force 
condition and the definition ([55]) of the excess potential 
(i.e., potential corrected for the buoyancy contribution) 
imply that the excess force and total force acting on the 
colloidal particles are directly related, 



V s * c = (1 - 0c) V s * c 



(69) 



Thus, for a uniform system, Eq. (|57p can be rewritten as 
the mobility relation 



where 



(70) 



(71) 



is the mobility coefficient describing the response of the 
system to the total (rather than excess) force. For com- 
patibility with the usual definition of the collective mo- 
bility in bulk suspensions, we present our results for the 
coefficient v instead of v 1 . 



Reflection method 



where k± is independent of z, owing to the continuity 
of the transverse momentum flux. Integrating the above 
relation with respect to z, and recalling that the effective 
viscous stress is associated with the excess part (|A13I) of 
the stress tensor, we find 

k = h(R\\ — Kj_ ) . (65) 

In our numerical calculations, the deviatoric part of the 
average stress tensor (|64|) is determined directly from the 
stresslet component of the induced-force distribution on 
the particle surfaces. To obtain the relation that links k 
to our numerical results, we rewrite equation l|64p. aver- 
aged across the film, in the form 

r' = a\+f3E, (66) 

where 

ol= |(kx + 2/C||), j3 = |(k|| - k±), (67) 



To obtain the shear and expansion film-viscosity co- 
efficients fj and k, the stress response of the film with 
a constant thickness h to the flow fields (|60p and l[6"3"[) 
has to be evaluated. Similarly, the collective mobility 
coefficient v is obtained by calculating the particle flux 
produced by a constant lateral force acting on the par- 
ticles. To solve these hydrodynamic problems we use 
a standard Stokcsian-dynamics algorithm for free space, 
combined with a flow reflection method. 

Our reflection technique is based on the well-known 
observation that the flow reflected from a free interface 
at z = has the form 

v out (x, y, z) = (l s - e z e z ) ■ v in (x, y, -z), (72) 

where v ln is the incoming flow. Accordingly, the image 
of a spherical particle centered at (X, Y, Z) and mov- 
ing with the velocity (U x , U y ,U z ) is a sphere centered at 
(X,Y,—Z), moving with the velocity (U x ,U y , — U z ). This 
reflection geometry is illustrated in Fig. [51 a). 
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To obtain the solution for particle motion in a film with 
periodic boundary conditions in the lateral directions, we 
use a Stokesian-dynamics algorithm for a periodic system 
of particles in free space. We construct an elementary cell 
in 3D from an elementary cell in the film and the reflec- 
tion of this cell in the upper interface, as illustrated in 
Fig. Multiple reflections of the composed elemen- 

tary cell in the 3D system ensure that the free-boundary 
boundary conditions are satisfied both on the lower and 
the upper film interface. Thus, the periodic 3D system 
and the doubly periodic liquid film with free interfaces 
are hydrodynamically equivalent. 

We note, however, that our hydrodynamic calculations 
in the thin-film geometry differ from the usual results 
for 3D periodic systems. Namely, in the thin-film case, 
the spheres are not allowed to overlap the lateral planes 
separating periodic cells, because these planes represent 
the film interfaces. This constraint is not required in 
simulations of suspensions in an unbounded space. 

In our numerical calculations we use the 
hydrodynamic-interactions algorithm introduced in 
[23| . with further improvements described in [24[. 
The algorithm combines a multipolar expansion of 
the flow field in the system, with a superposition of 
two-particle lubrication corrections. Our algorithm 
includes the multipolar expansion of Stokes flow to an 
arbitrary order, thus accurate results, with controlled 
multipolar-truncation error, can be obtained. 

C. Simulation details 

The short-time film viscosity coefficients f) and k were 
determined using relation 

N 

f' = V- l Y,Si (73) 

between the average stress tensor ([52"]) or (|rJoT) and the 
stresslets 

Si = ~J [r 2 F 4 (r,) + F^rOr* - fr, • F 4 (r. ( )l] dr t (74) 

of the induced-force distribution Fj on the surfaces of 
particles i — 1, . . . , N (where Yi = r — Ri is the rela- 
tive position vector with respect to the center of particle 
i) . The film- viscosity coefficients f] and k were evaluated 
using relations (|62|) and ([68]) for a suspension of spheres 
in the external 2D straining flow (|60|) and axisymmetric 
straining (|63|) flow, respectively. The collective lateral 
particle mobility was obtained from the relation 

N 

v = \ N ~ l Z>« :l »' ( 75 ) 

where /i" is the translational component of the many- 
particle translational mobility matrix in a periodic sys- 
tem with zero net flow. 




h/d 

FIG. 7: Normalized particle contributions to the (a) shear 
and (b) expansion film- viscosity coefficients, as defined by 
Eqs. (|76[) . versus film thickness h, for particle volume frac- 
tions 4> c — 0.1 (dotted line), 0.2 (dashed), 0.3 (dash-doted), 
0.4 (solid). 

The quantities (|73|) — (|75|) were determined for an equi- 
librium distribution of N = 40, 80 and 160 spheres in 
films of different thickness h. The results were averaged 
over several hundred particle configurations, obtained us- 
ing a Monte-Carlo method. In our calculations, the force 
multipoles with the vector-spherical-harmonic orders up 
to I = 5 were included. 



D. Numerical results 

Our numerical results for the shear and expansion film- 
viscosity coefficients fj and k are depicted in Fig. and 
those for the collective lateral mobility coefficient v are 
shown in Fig. [8] The film-viscosity coefficients fj and k 
are presented in terms of the normalized particle contri- 
butions f] p and k p , defined by the relations 

V = Voh(l + §0 C 77 P ), (76a) 
k = ^-<f> c T)QhR p , (76b) 
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FIG. 8: Normalized collective mobility coefficient, as defined 
by Eq. (|77[) . versus film thickness h, for the same values of 
volume fraction <f> c as in Fig. [7] 



where 770 is the viscosity of the suspending fluid. Particle 
mobility in the film is characterized by the dimensionless 
mobility coefficient v v , normalized by the mobility of an 
isolated particle in free space 



v = STtrjodiSp 



(77) 



With the above definitions, we have fj p , k p , v p — > 1 in 
the low-density and thick-film limit <p c — > and h/d —> 
00. The transport coefficients fj p , k p and v p are plotted 
in Figs. [7] and [5] versus the film thickness h for several 
particle volume fractions 4> c . 

The results in Fig. [7] indicate that the normalized film- 
viscosity coefficients depend on the film thickness in an 
oscillatory way, especially in the high volume-fraction 
regime. The oscillations are more pronounced for the 
shear- viscosity coefficient, which has a large peak at 
h/d w 1.4 [close to the position of the peak of the lateral 
pressure p c ii, shown in Fig. [5^6)]. At high volume frac- 
tions, the normalized mobility coefficient v p , depicted in 
Fig. [8j also has significant oscillations. The oscillations 
are more clearly visible in the plot shown in Fig.[S]Ja) on 
an expanded scale. 

According to the results depicted in Fig.JTJa), the nor- 
malized shear viscosity of the film f) p assumes a minimal 
value at h/d = 1. This behavior stems from the fact that 
for lateral particle motion, the image particles move in 
the same direction as the original ones (cf. the schematic 
representation in Fig. [6J. Such a collective motion results 
in a reduced dissipation. For the same reason, the col- 
lective mobility v p has a maximum at h/d = 1, as shown 
in Fig. [51 In contrast, the expansion- viscosity coefficient 
tip, plotted in Fig. 0J&), has an 0(e _1 ) singularity at 
h/d = 1 (where e = h/d — 1 denotes a typical gap be- 
tween the particle surfaces and the film interfaces), owing 
to the lubrication resistance between the original parti- 
cles and their images (which have the opposite normal 
velocity component). At low particle volume fractions 
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FIG. 9: Normalized (a) collective mobility coefficient and (6) 
particle contribution to the expansion viscosity of a particle- 
stabilized film, versus film thickness h, for <f> c = 0.4 and 
spherical-harmonics truncation order Z max = 2 (dotted line), 
3 (dashed) 4 (dash-dotted), and 5 (solid). 



the collective mobility coefficient v p diverges as log(</j c ), 
because of the logarithmic far-field behavior of 2D Stokes 
flow. 



E. Convergence tests 

The results for the film-viscosity and mobility coeffi- 
cients shown in Figs. [7] and [8j were obtained for a system 
of TV = 80 particles using induced-force multipoles up 
to the spherical-harmonics order Z max = 5. The conver- 
gence of the results with ^ max is illustrated in Fig. [5] and 
with the system size N in Fig. [10] for the highest particle 
volume fraction considered in our simulations, <f> c = 0.4. 

The results shown in Fig. [D^a) for the mobility coef- 
ficient v indicate that the truncation of the multipolar- 
expansion at the order l max = 2 yields an error of 20 % to 
30 % for this quantity. The results obtained with Z max = 3 
have the accuracy better than 8% and Z max = 4 yields 
results accurate within 1 %. We note that Z max = 2 in- 
cludes all multipoles of the STD scheme used by Brady 
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FIG. 10: Normalized particle contributions to the shear vis- 
cosity of a particle-stabilized film, versus the film thickness h, 
for 0c = 0.4 and system size N — 20 particles (dotted line), 
40 (dashed), 80 (dashed-dotted), and 160 (solid). 

and his collaborators [3l[ , and thirteen other flow multi- 
poles (including potential flows). Thus, for the mobility 
coefficient v the standard STD scheme is insufficient for 
obtaining accurate results, even if the lubrication correc- 
tions are included in the numerical scheme. 

Our numerical results for the film- viscosity coefficients 
f) and k indicate that the convergence with Z max is faster 
for these quantities than the convergence for the collec- 
tive mobility v. As illustrated in Fig.[9^&), the truncation 
at Z max = 2 already yields quite accurate results for the 
coefficient k, except for films with h/d < 1.5. For fj the 
convergence with Z max is similar. 

The convergence of our results with the system size is 
illustrated in Fig. \W[ where the results for the effective 
shear-viscosity coefficient fj, evaluated with l max — 2, are 
shown for systems with N between 20 and 160 particles 
in a unit cell. The results indicate that the convergence 
of the shear viscosity fj is fast. A similar convergence has 
been obtained obtained for the expansion viscosity k and 
the collective mobility v. 



V. CONCLUSIONS 

Our unified quasi-two-dimensional description of equi- 
librium states and transport properties in particle- 
stabilized thin liquid films relics on the film-tension rep- 
resentation of the film thermodynamics. In this repre- 
sentation the system is characterized by extensive pa- 
rameters corresponding to conserved quantities: the film 
volume, film area, and the number of fluid molecules and 
colloidal particles. The conjugate intensive parameters 
are the normal pressure, film tension, and chemical po- 
tentials of the fluid phase and the colloidal component. 

Equilibrium conditions and their physical interpreta- 
tion have been analyzed for nonuniform films in external 



potentials. In particular, we have shown that the particle 
contribution to film tension is associated with the excess 
lateral pressure beyond the uniform isotropic background 
pressure p( c \ Since the normal pressure in the film is 
equal to the external pressure, the excess lateral pressure 
equals the difference between the lateral and normal pres- 
sure components. 

Based on our thermodynamic formalism and the local- 
equilibrium assumption, we have also introduced a quasi- 
two-dimensional description of transport processes in 
particle-stabilized liquid films. In some aspects, a film 
composed of an incompressible suspension is analogous 
to a 2D compressible fluid, with the film thickness play- 
ing the role of the mass density per unit area, and film 
tension being similar to the pressure (with a minus sign). 

When deviations from equilibrium are small, the mo- 
mentum transport in the film is characterized by two lin- 
ear viscosity coefficients: the shear and expansion viscosi- 
ties. The particle diffusive transport is described by two 
kinetic coefficients relating the particle flux to the gra- 
dients of the normal component of the osmotic pressure 
p c ± and the excess particle chemical potential /i' c . Using 
the Gibbs-Duhem relation in the film-tension representa- 
tion, the two kinetic coefficients have been combined into 
a single collective mobility coefficient. Our explicit calcu- 
lations of the short-time linear transport coefficients for a 
film stabilized by a hard-sphere suspension indicate that 
these coefficients are oscillatory functions of film thick- 
ness. The oscillations are due to the short-range layering 
order in the film. 

While the results of our investigations provide a conve- 
nient framework for describing equilibrium and nonequi- 
librium properties of particle-stabilized thin liquid films, 
there still remain many open questions. First, we consid- 
ered only films with surfactant-free interfaces. Surfactant 
adsorbed on the interfaces would contribute to film ten- 
sion, but otherwise it should not significantly alter the 
equilibrium conditions in the film. It would, however, 
substantially influence the relaxation times for decay of 
nonuniformities of different intensive thermodynamic pa- 
rameters. Thus partial equilibrium states in stratified 
films (discussed in [l7|) can be significantly affected. 

Adsorbed surfactant modifies not only the equilibrium 
equations of state, but also changes the film viscosity 
and particle mobility, because of its effect on the bound- 
ary conditions at the film interfaces (c.f., our analysis 
of this problem presented in [32|). The hydrodynamic- 
interaction algorithm used in the present paper is ap- 
plicable only to surfactant-free interfaces. However, 
our recently developed Cartesian-representation method 
[HlHl, HI]) can be applied to evaluate the effective trans- 
port coefficients in a particle-stabilized film with more 
complex boundary conditions. 

Another important problem that has not been ad- 
dressed here is transport through contact lines between 
phases of different thickness in stratified films. An anal- 
ysis in [36[ suggests that such transport may control 
the rate of expansion of thin spots in the film during a 
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stepwise-thinning process. On the other hand, hydrody- 
namic scaling arguments [l7| seem to imply that particle 
transport through a contact line should be faster than 
particle diffusion (at least in surfactant-free films). This 
transport, however, can be significantly slowed down, ow- 
ing to surfactant effects, and because of the entropic bar- 
rier associated with a change of particle arrangement at 
the line separating phases of different thickness. This 
problem will be analyzed in our future publications. 
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APPENDIX A: EFFECTIVE SURFACE STRESS 

In this appendix we derive the effective lateral stress- 
balance equation under creeping-flow conditions in a thin 
film with varying width h. The film occupies the region 
between the lower and upper interfaces 

z = z L (x,y), z = z v (x,y). (Al) 

The region occupied by the film can be described by the 
characteristic function, 

Of (r) = 6{z - z L )6( Zl j - z) (A2) 

where 9 denotes the Heaviside step function. 

The film is surrounded by a gas with a constant pres- 
sure p^ , and a lateral force field 

f(r) = f (r)f(r) (A3) 

acts on the film. Accordingly, the stress tensor in the 
whole 3D space 

T t ot = -(l-0f> (c) l + 0fT (A4) 

(where r denotes the stress tensor in the film, and p^ is 
the external pressure) satisfies the full 3D balance equa- 
tion 

V • Ttot = f • (A5) 

There are no capillary-pressure and Marangoni-stress 
contributions in the above equation, because of our as- 
sumptions that the curvature of the interfaces is small, 
and the interfacial tension is constant. If needed, such 
contributions can be added as appropriate delta distri- 
butions of the force on the film interfaces. 



To obtain the film-tension representation for the lateral 
stress balance in the film, it is convenient to reorder terms 
in Eq. (fM)) . 

xtot = -p^l + BtSr. (A6) 
The first term on the right-hand side of this equation 
is the constant background pressure that does not con- 
tribute to the stress-balance (|A5[) . The second term rep- 
resents the excess stress in the film. Using the continu- 
ity of the transverse momentum flux and the associated 
boundary condition (|5ip for the normal stress we find 

5t = (t|| - t ± )\ s . (A7) 

By inserting (|A6[) into the stress-balance equation (|A5|) 
and integrating with respect to the variable z, we obtain 

J dz (e z ^- + V.) • 6 f Sr = -f, (A8) 

where f is the force density per unit area of the film. The 
first term on the left-hand side of the above relation van- 
ishes because 0{(x,y,z) = for z = ±oo. In the second 
term we change the order of integration and differenti- 
ation to obtain the quasi-two-dimensional stress-balance 
equation for the excess stress, 

V s • hSt = -f, (A9) 

where the bar denotes the average $2$. 

In an equilibrium state, only the equilibrium pressure 
tensor (fTJ contributes to t, 

t = -(p±e z e z +P||I„). (A10) 

Taking into account the boundary condition ([7]) and the 
definition of film tension ([5]) , we find that relations (|A9[) 
and (|A10[) yield the lateral stress balance (fTT|) . 

In a more general, nonequilibrium case, Eq. (|A9|) is 
split into the lateral and vertical components, 

V s -Mf-l s = -f. (All) 

V s -hST-e z = 0, (A12) 

Since the normal component (|A12|) vanishes on the as- 
sumption that there is no transverse external force, we 
find that (|A11I) yields the lateral stress balance (|4T|) . with 
the effective 2D stress 

T = h5f-\ s . (A13) 
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